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Abstract 

The thermodynamical stability of DNA minicircles is investigated by means of path integral 
techniques. Hydrogen bonds between base pairs on complementary strands can be broken by 
thermal fluctuations and temporary fluctuational openings along the double helix are essential 
to biological functions such as transcription and replication of the genetic information. Helix 
unwinding and bubble formation patterns are computed in circular sequences with variable radius 
in order to analyze the interplay between molecule size and appearance of helical disruptions. The 
latter are found in minicircles with < 100 base pairs and appear as a strategy to soften the stress 
due to the bending and torsion of the helix. 
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It is known that the helicoidal conformation of DNA is essentially determined by the 
hydrophobicity of purine and pyrimidine bases and by the bond angles in the flexible sugar- 
phosphate backbone while sequence of the bases and environmental conditions due to the 
solvent also contribute to the molecule shape [1]. As each strand bears a negative charge 
(e) for each phosphate group and the rise distance between adjacent nucleotides is rs./ 3.4A, 
the bare double helix has a high linear charge density of ~ 0.6 e A -1 . Although the effective 
charge density is reduced by the counterions in the solvent, the electrostatic strands repulsion 
is key to the stability of the helix and also affects the inter-helical chiral interactions in those 
condensed phases of DNA assemblies (such as liquid crystals) which underlie the impressive 
growth of DNA-based structures recently witnessed in materials science. 

Base pairing and base stacking are the fundamental interactions which control the syn¬ 
thesis of DNA and determine the thermodynamic stability both of single helices and of 
helix aggregates. However even stable duplexes at room temperature show local openings, 
temporary bubbles, which are intrinsic to the biological functioning as they permit the 
transcription and replication of the genetic code. Such bubbles are due to the strong fluctu- 
ational effects on the hydrogen bonds between complementary strands and cause the local 
unwinding of the double helix which ultimately leads to a state of negative supercoilcd DNA 
for almost all living beings. While these processes are qualitatively understood, quantitative 
predictions of (energetically) optimal helical configurations for specific systems are scarce. 
We contribute to fill this gap by introducing a new path integral computational method 
which readily applies to loops of DNA as those found in bacterial plasmids, viral genomes 
and also mammalian cells. 

I. MESOSCOPIC MODEL 

Let’s arrange N base pairs (bps) on a circle with radius R such that, 2ixR/N ~ 3.4 A, as 
depicted in Fig. [T| When all bps centers of mass lie on the circumference, which represents 
the molecule backbone, the system is in the ground state. 

Say q, the inter-strand fluctuation for the Abase pair (i = 1 ,..N) with respect to the 
ground state. Hence, we define the vector tp 

(ti) x = |q| cos fa cos Of, (U) y = (R + |q| sin#*) cos fa] (ti) z = (R + |r;|) sin fa . (1) 
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The ground state is recovered once all ^-fluctuations vanish hence, U = R, Vi. The 
polar angle, 9i = (i — 1)27 r/h + 9s, measures the i— bp twisting around the molecule 
backbone, with h — N/Tw being the number of bps per helix turn and Tw is the twist 
number accounting for the coiling of the individual strands around the helical axis l|. The 
azimuthal angle, (pi = (■i —1)27 r/N + (ps, defines the bending between adjacent bps along the 
stack. As the polynucleotide chain has a direction due to the chemistry of the intra-strand 
bonds, a distribution of values for the twist (6s) and the bending (cps) of the first bp in 
the sequence is weighed in the computation. The fluctuational orbits defined by i — 1 and 
i = N + 1 overlap consistently with the closure condition holding for the DNA ring. 


II. SPACE-TIME MAPPING 


My previous path integral analysis of DNA \2\ ll()| were based on the ansatz that the 
bps displacements could be treated as one dimensional paths x(tP), || —* xfa), with the 
imaginary time r* G [0, /3] and (3 being the inverse the temperature. 

Here I introduce a more general (albeit more CPU time consuming) space-time mapping 
technique, which does not pin a base pair to a specific r* thus avoiding the somewhat 
arbitrary partition of the f3 length in N intervals: 


|u| —>■ Xi(r); |U| -A ^(t); tG[0,/?]. (2) 

The paths Xi(r) are expanded in Fourier series, the constraint ay=i(r) = Xi= n+i( t ) f° r 
a DNA loop is implemented and, for the model in Fig. |TJ the bps fluctuations in the path 
integral formalism are given by: 


V iij) = [R 2 + Xi(r) 2 + 2R\xi(T)\f(6i,(pi)] 1/2 

f(9i, (pi) = sin 9i cos 2 (pi + sin 2 <pi. (3) 

This new mapping technique lets r as an integration variable and permits to set a re¬ 
alistic rise distance between adjacent nucleotides along the stack. Written in terms of the 
r/,(r), the Hamiltonian contains: i) a Morse potential, Vm\j]i(t)\, describing the effective 
hydrogen bond interaction between bps mates, ii) a solvent term, V so i [r/,;(r)], accounting for 
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FIG. 1: (Color online) (a) Helicoidal model for circular DNA with bending planes. The blue filled 
circles are the centers of mass of the base pairs stacked along the molecule backbone with a rise 
distance 3.4 A. In the ground state all bps lie on the circumference with radius R. The red-shaded 
areas are spanned by the fluctuational vectors whose amplitude is measured by |?’i| for the i— bp. 
The azimuthal angle (j>i measures the bending of the i— bp plane with respect to the (x 1 , y') plane, 
x’ being normal to the sheet plane, (b) Local reference system for the i— bp. 9i is the twist around 
the molecule backbone. The z-axis is tangent to the ground state circle. 


hydrogen bond recombination with the counterions dissolved in water, in) a (two particles) 
stacking potential, Vs [rji (r), 7^-1 (t)] , between adjacent bases along the strand. The poten¬ 
tial parameters are widely discussed in Refs. jd-6]. Accordingly, the classical action A[ry] 
is a T-integral of the base pair kinetic plus potential energies [9| and the partition function 
Z N for a sequence with N bps reads: 


4 














( 4 ) 


N „ 

zn= n f Dxi exp { 

i=l 9s,4>s 

N rP 

A fo) = 

i= 1 7 ° 


%Vi( T ) + V M [Vi(T)\ + V ao i[r}i{r)] + V s [^(t), Vi-i^)] 


where /i = 300 amu is the bp reduced mass. The measure Dxi is a multiple integral over 
the path Fourier coefficients while § indicates that the paths Xj(r) are closed trajectories. 
The integration over the two particles potential greatly enhances the computational times 


with respect to the previous method 
stacking in sequences of any length. 


TO] but it offers a more realistic model for the base 


III. FREE ENERGY OF MINICIRCLES 


There has been a widespread interest in the properties of small sequences of DNA fol¬ 
lowing measurements of high cyclization probabilities which have pointed to an intrinsic 
flexibility for fragments of ~ 100 bps or less ll|. It has been suggested that the formation 
of local disruptions in the helix of DNA loops may be the mechanism which permits to re¬ 
lease the torsional stress and energetically favors the stability of bent molecules. Certainly, 
the breaking of some bps and the opening of fluctuational bubbles change the average pitch 
of the helix, that is the number of bps per helix turn. 

In the theoretical framework synthesized by Eq. (J4]), I compute the free energy (F = 
/3^ 1 \uZn) of three heterogeneous loops with different length, N = 18, 66, 86, but similar 
content of AT -bps and GC -bps, ~ 50% each. The rise distance is pinned to the experimental 
value. The 66— and 86— bps loops have been prepared as described in Ref. [12] whereas 
the N — 18 bps is a toy sequence here introduced for comparison. For each loop, we 
simulate a broad set of values for the twist Tw and, for any Tw (i.e., h), Z N is calculated by 
summing over an ensemble of fluctuational paths representing a large number of molecule 
configurations, about 10 7 paths for each base pair in Eq. (Iffi). 

Hence, the free energy is obtained as a function of the helical repeat h. The computational 
time for a simulation, e.g. for the N — 66 sequence, is about eight days on a workstation 
(Intel Xeon E5-1620 v2, 3.7GHz processor). The room temperature results for F/N are 
plotted in Fig. [2j While the shortest loop shows a free energy minimum also at h ~ 7, the 
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FIG. 2: (Color online) Free energy per base pair for three loops with N = 18, 66, 86. The free 
energies are computed, at room temperature, as a function of the number of bps per helix turn. 


minima for all loops are remarkably found for h in the range ~ (10 —12) in line with the well 
known values of helical pitch in DNA sequences jl]. Even more interestingly, the 66— and 
86— minicircles also show free energy minima at larger h due to a spontaneous unwinding of 
the complementary strands. This Tw reduction is consistent with the observations 12] of bps 
disruptions occuring in small loops as a consequence of the strong rotational deformations 
of the helix. The bending stress decreases by increasing the radius of the loop. 

Then, for an ensemble of molecules with N bps, the free energy minimization evaluates 
the thermodynamically stablest values of helical repeat, that is an average property of the 
molecule ensemble. Furthermore, the path integral method can also determine the proba¬ 
bilities for the formation of fluctuational bubbles and select those base pairs along the DNA 
sequence for which hydrogen bond breaking is more likely to occur. By tuning the system 
temperature and the input parameters which control the counterions concentration in the 
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solvent, we thus obtain a general and reliable computational scheme for the modeling of 
heterogeneous DNA loops in various ambient conditions. 
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